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Abstract. - We show how the conversion of magnetic field energy via magnetic reconnection 
can progress in a fully three-dimensional, fast, volume-filling regime. An initial configuration 
representative of many laboratory, space and astrophysical plasmas spontaneously evolves from 
the well-known regime of slow, resistive reconnection to a new regime that allows to explain the 
rates of energy transfer observed in jets emitted from accretion disks, in stellar/solar fiare processes 
as well as in laboratory plasmas. This process does not require any pre-existing turbulence seed 
which often is not observed in the host systems prior to the onset of the energy conversion. The 
dynamics critically depends on the interplay of perturbations developing along the magnetic field 
lines and across them, a process possible only in three-dimensions. The simulations presented 
here are the first able to show this transition in a fully three-dimensional configuration. 



Introduction. — A continuous and volume-filling pro- 
cess converting magnetic energy into high energy parti- 
cles is manifestly present in jets emitted from accretion 
disks around forming stars or black holes. In fact the 
synchrotron emission from the whole jet means that an 
active process is actually sustaining and reaccelerating a 
population of relativistic particles that hence emit. In 
absence of such a mechanism the relativistic population 
would be confined close to the energy source (e.g., the ac- 
cretion disk) and quickly lose its energy by radiation as 
it moves away. The sustained production and accelera- 
tion of extremely energetic particles in vast macroscopic 
regions is also observed in solar system plasmas (Sun and 
Earth environments) and in laboratory. The mechanism 
underlying the above mentioned phenomena is yet to be 
fully understood, but it must clearly include two critical 
features: the acceleration process must be volume-filling 
and it has to be sufSciently fast that the energy input from 
this mechanism matches the observed rates of energy con- 
version in Nature. 

Understanding the process of fast magnetic energy con- 
version remains a key task in plasma physics. Magnetic 
reconnection is the leading process long suggested to ac- 
complish this task. In magnetic reconnection, field lines 
break and change their topological connection in a way 



that makes a piece of a field line connect with a piece 
of another field line [lj[2]. In the process the magnetic 
energy is released in the form of kinetic energy, that is 
acceleration and heating. However, under any understood 
mechanism, on macroscopic scales, reconnection requires 
a dissipation process such as resistivity allowing the con- 
version, but in all observed natural and in most man-made 
conditions the level of dissipation measured or induced is 
several orders of magnitude too low. Also so-called anoma- 
lous dissipation processes have been traditionally invoked. 
However, the microscopic mechanisms leading to the for- 
mation of the invoked anomalous effects depend on several 
factors and have not yet been conclusively identified. For 
example within the solar corona, the free non-potential en- 
ergy of the (local) solar magnetic field is released via mag- 
netic reconnection both during sudden and violent erup- 
tive processes such as fiares and coronal mass ejections 
(CME) or it can occur in a slowly quasi-steady way which 
may contribute to the coronal heating. Yet measured high 
temperatures of current-sheets formed in the aftermath of 
CME's [sjjl] and inferred characteristic dynamic (Alfven) 
time scales throughout the different phases of solar ex- 
plosive phenomena have pointed out a number of discrep- 
ancies between the level of plasma resistivity required in 
the solar corona, e.g. in the case of stationary current- 
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sheets [5] wherein it is about 10^^ — 10^^ m^s"^, and the 
one deduced from theories of classical and anomalous re- 
sistivity, that is about 1 m^s~^ and 10^ — 10^ m^s""'^, 
respectively. This implies that a very efficient diffusion 
mechanism is at work and that neither the classical resis- 
tivity nor anomalous effects are necessarily governing this 
process. 

At the kinetic level, which is measured quantitatively 
as the radius of gyration of an ion particle in the ambient 
magnetic field, reconnection can be locally fast without re- 
quiring the dissipation mechanisms needed on macroscopic 
scales. These scales have been recently analyzed especially 
for the jet originating from the relatively distant radio 
galaxy 3C303 (z — 0.141) [6|. The jet presents a knotty 
structure, that is brightly radiating spots alternate with 
dim segments, observed in both radio and X-ray bands and 
these knots may be considered as actual plasma features 
created near the central source and propagating along the 
jet axis. In this case the mechanism at work could be the 
reconnection of magnetic field lines [jj. It is conceivable 
that due to some yet unidentified mechanism kinetic-scale 
reconnection might be allowed to fill the whole jet and pro- 
vide energy conversion at the macroscopic scales needed. 
However, very recent works have suggested that it is not 
necessary to invoke such extreme extrapolation from the 
observed scales of energy conversion to the minute scales 
where kinetic reconnection becomes active. We show here 
the first demonstration via numerical simulations of the 
ability of macroscopic processes of magnetic reconnection 
to achieve the two required goals needed to explain the 
observations: the process is fast enough and it is volume 
filling. We remark that both results have long eluded the 
research. Fast reconnection within macroscopic models 
described by the viscous-resistive magnetohydrodynam- 
ics (MHD), which does not include kinetic processes or 
anomalous effects, have long been believed impossible. 
The consensus was that under these conditions reconnec- 
tion would progress under the so-called Sweet-Parker (SP) 
rate and hence determine a very slow reconnection. This 
point of view changed in the last couple of years due to 
the discovery of a new regime where reconnection is un- 
steady and self-driven [8 -10 . The typical state involved 
in the SP mechanism is destabilized and a new turbulent 
layer forms in its place. Yet due to the limitation of be- 
ing two dimensional these discoveries still failed to explain 
the volume- filling feature of reconnection and according to 
them the fast reconnection layer remained localized. 

A possible way to overcome this last step was indi- 
cated by Lazarian and Vishniak who proposed that a 
fully three dimensional reconnection layer where turbu- 
lence is included could provide a macroscopic and fast 
reconnection scenario 11 . This idea had also been pro- 
posed in a cartoon-based speculative scenario by Shibata: 
the current-sheet undergoes recurrent secondary tearing 
instabilities 12 ■ 15 occurring at different spatial scales 
and this leads to a fractal (i.e., scale- free) current-sheet 
structure, where many magnetic islands of different sizes 



(subject also to coalesce with each other) connect macro- 
scopic and microscopic scales and fill the scale gap men- 
tioned above 



16 



The main idea is that the energy is 
supplied at the largest spatial scales and it is being re- 
distributed via MHD turbulent cascade to smaller scales 
and finally dissipated. The first attempts to study this 
process were done in two dimensional simulations on the 
effect of MHD turbulence on magnetic-reconnection 
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and nowadays not only it is suggested to be a fundamen- 
tal mechanism to model solar explosive phenomena 18 



but also recently it has been claimed to have obtained the 
first in situ evidences of magnetic reconnection in turbu- 
lent plasmas within the heliospheric environment 
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In the present paper we demonstrate quantitatively that 
fast, three dimensional and volume filling reconnection 
can be achieved spontaneously within the framework of 
the pure viscous-resistive MHD and even without invok- 
ing any externally driven turbulence. In fact the simu- 
lations presented here have an unprecedented low level 
of numerical diffusion that in previous studies eliminated 
the transition to the fast turbulent regime (replacing it 
with a fast laminar regime purely driven by numerical ef- 
fects). Indeed a system completely laminar in conditions 
similar to those found in many astrophysical and labo- 
ratory systems is observed to settle at first in the well 
know slow SP regime, but its subsequent evolution spon- 
taneously brings it through three stages. First, the desta- 
bilization of the SP regime previously studied in two di- 
mensions repeats itself and even in three dimensionsional 
simulations the system still maintains a two-dimensional 
behavior. Secondly, an explosive growth of instabilities 
in the third dimension breaks the translational invariance 
and makes the evolution fully three-dimensional. In the 
third stage, a new regime is eventually observed, where 
the reconnected states are fully mixed and the process 
encompasses the whole domain. The transition in 3D, 
unlike in 2D, allows to reach regions of macroscopic size 
and providing the needed macroscopic channel to explain 
the observed amounts of energy conversion. Furthermore 
the dynamics critically depends on the interplay of per- 
turbations developing along the magnetic field lines and 
across them, a process possible only in three-dimensions. 
Both the volume-filling and the fast-progressing features 
of magnetic reconnection are within the same model and 
without requiring extra physical processes such as anoma- 
lous resistivity. 

Method. — In order to reach the needed high resolu- 
tion and to reduce the influence of spurious numerical ef- 
fects, we use a code with a fully-implicit (in time) particle 
(or FLIP) algorithm to solve an Eulerian-Lagrangian for- 
mulation of three dimensional resistive and viscous MHD 
equations, FLIP3D-MHD 



20 . This algorithm solves the 



MHD equations in four steps and it uses both a particle 
and a grid description of the fluid: 

1. particle data are interpolated to a Lagrangian grid 
co-moving with the fluid; 
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2. Lagrangian finite difference equations are solved on 
the moving grid to advance the fields; 

3. changes in the grid data over a time step are interpo- 
lated to the particles to advance the particle data; 

4. a new grid is generated (to avoid the tangling of a 
Lagrangian mesh) and the particles are interpolated 
to the grid for the next time step 

The most important difference between FLIP and earlier 
particle-in-cell methods 21 is in step 3. By interpolat- 



ing changes in the solution from the grid to the particles, 
rather than by replacing the particle data with the new 
grid solution. Consequently, computational diffusion de- 
creases with the time step. At. This change guarantees 
consistency, even with higher order interpolation, and re- 
duces computational diffusion below that achievable with 
Eulerian finite difference methods with a similar range of 
applicability. 

The particle-in-cell (PIC) method is Galilean invari- 
ant 22 and FLIP extends this property to the MHD equa- 



tions. Galilean invariance is achieved by making the mag- 
netic field a particle variable. Changes in the magnetiza- 
tion are interpolated to the particles to calculate the new 
particle magnetic moment. This last step introduces some 
computational diffusion. However, the diffusion is small 
compared with that produced by, e.g., van Leer convec- 
Spatial differencing insures that differencing at 



20 



tion 

grid vertices and cell centers are conjugate operations [23 



Thus momentum is conserved, and energy is conserved 
with implicit differencing in time. 

FLIP has been validated by comparisons of analyti- 
cally solvable problems, such as the Rayleigh- Taylor in- 
stability, confined eddy, and Kelvin-Helmholtz instabil- 
ity 20 24 25 , and by comparisons with other numerical 



computations. 

In the simulations presented here, the visco-resistive 
MHD equations have a uniform viscosity and resistivity 
parametrized by the global Lundquist number, 5, and ki- 
netic Reynolds number, Rm, both set to 10^. We define a 
stream- wise direction {X) with periodic boundary condi- 
tions at X = and Lx, a cross-stream (Z) and span- wise 
(Y) direction with, respectively, reflecting and periodic 
boundary conditions at Z{Y) = and Z{Y) = Lz{Ly)- 
We consider 360 (in X) x 60 (in Y) x 240 (in Z) La- 
grangian markers arrayed initially in a 3 x 3 x 3 uniform 
formation in each of the 120 x 20 x 80 cells of our numer- 
ical grid. The mathematical formulation of the problem 
include an ideal equation of state with polytropic index 7 
such that p = p(7 — 1)/ where p, p and I are respectively 
the kinetic pressure, the density and the enthalpy of the 
plasma. The basic formula has the main field component 
reversed, Bx{Z) — tanh(Z/L), over a characteristic 
length scale L. Additional component can be considered 
to obtain the so-called guide field configuration (with a 
uniform magnetic field component in the direction orthog- 
onal to the gradients and the main field), By{Z) = aSo 
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Fig. 1: General evolution of the current-sheet dynam- 
ics: electric field and different contributions to the en- 
ergy during the simulation performed with FLIP3D- 
MHD. Panel (a) shows the parallel electric field, ^par ~ 
E • B/|B|. Panel (b) shows the different magnetic and kinetic 
contributions to the system energy during current-sheet recon- 
nection dynamics. Moreover the system is sampled at six re- 
markable instants indicated by ti = 72, t2 = 96, the time 
interval between ts = 104 and t4 = 112, ts — 128, and t = 160 
clearly shown in the figure. 



(a — constant), or the force-free configuration characteris- 
tic of situations where the plasma is at its minimum energy 
state and all forces are absent, By{Z) = i?o/cosh(Z/L). 
While all cases were considered in the present study, the 
results will focus on the simples case of By = 0. 

Sv^reet-Parker layer and its destabilization: mag- 
netic connectivity and energy evolution. When 
regions of different magnetic field orientation come in con- 
tact, the simplest case being two regions of opposite polar- 
ity, the possibility arises for the process of reconnection. 
For astrophysical jets, this could be the case of the field 
generated from the accretion disk and carried along the jet 
by the Poynting flux. Magnetic field lines coming from the 
disk must return to it (assuming the classical absence of 
magnetic monopoles) , leading to regions of different mag- 
netic field orientation. Moreover the rotation of the disk 
transmits a rotation also to the field lines, which are hence 
wound around the jet, and this facilitates the insurgence 
of regions of reconnection 



26 . Similar situations are ob- 



served in the solar system as well, i.e. regions of opposing 
polarity coming in contact in the photosphere of the Sun, 
opposing magnetic field lines in the solar wind or in the 
magnetosphere (during periods of southward directed in- 
terplanetary magnetic field compared with the northward 
directed Earth field) , or the Earth field lines bent and push 
together in the Earth magnetic tail. In laboratory, such 
configurations are also present by design in fusion devices 
and in devices especially designed to study the process of 
magnetic reconnection [27]. 

In the panel (a) of Fig. [T] the parallel electric field 
(= ^pAR = E • B/|B|) integrated in the whole numerical 
domain is shown as a function of time: the rate of changes 
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Fig. 2: Volume iso-surface of the density in the region 
of the current-sheet at the centre of our numerical do- 
main and magnetic field lines. The blue-red color table 
for the magnetic field lines describes the intensity of the X- 
component. The panel (a) shows the system at t ~ 72 (ti) 
when the disruption in multiple pieces of the current-sheet 
starts. The panel (b) shows the system at t ~ 96 (t2) when 
coalescence process leads to the formation of an unique plas- 
moid at the centre of the box. 




Fig. 3: Volume iso-surface of the density in the region 
of the current-sheet at the centre of our numerical do- 
main and magnetic field lines. The blue-red color table 
for the magnetic field lines describes the intensity of the X- 
component. The panel (a) shows the system at f « 104 (ts) 
and in the panel (b) a detailed clip of the central magnetically- 
confined higher-density region wherein It is evident the sinuous 
mode of the kinking dynamics of the structure and the com- 
plex closed loops of the magnetic field following the instability 
dynamics of the island. 



in the magnetic connectivity is associated with the maxi- 
mum value of $ PAR 28 . At t ~ 72 (dimensionless units of 



time measured with respect to the typical Alfvcn time) the 
system experiences an abrupt switch from an initial slow 
and laminar reconnection regime to a fast reconnection 
process. The parallel electric field evolution remarkably 
underlines the spontaneous transition to the fast regime 
lasting from t = 72 to about t — %Q where the (first) 
maximum of ^par is reached. The second, and highest, 
peak describes a third reconnection phase occurring now 
along the span- wise direction {Y), that is the direction 
orthogonal to the first tearing plane: it can be present 
only in fully 3D configurations and, as we show, it corre- 
sponds to a macroscopic turbulent configuration. In fact 
this last phase dispalys a volume-filling kinking instability 
that drives the system evolution towards a final chaotic 
state. In the panel (b) the evolution of the magnetic and 
kinetic contributions to the system energy is considered. 
The evolution is divided into several time intervals from 
ti = 72 to tg = 160 (dimensionless units of time) useful in 
the analysis of the subsequent data. The initial reservoir 
of magnetic energy, which corresponds essentially to the 
B\/2 contribution, is rapidly converted to kinetic energy. 

In Fig. [2][4] the volume iso-surface of the density to- 
gether with the magnetic field lines (colored according to 
the value of the X-component of the magnetic field) are 
shown for t « 72 (ti) and t = 96(^2) (see Fig. [l] for the 



definition of ti with i = 1,...,6). As shown in Fig. [2] 
panel (a), at t = ti, the current-sheet, thus far undergo- 
ing a slow Sweet-Parker (laminar) reconnection, destabi- 
lizes and tears into multiple-islands separated by smaller- 
scale current-sheets: the system chaotically evolves to- 
wards smaller and smaller length-scales, but it still main- 
tains a two-dimensional configuration on any X-Z plane of 
the numerical domain. Previous studies have analyzed this 
phase, its fractality, turbulence and its scaling [9j|29 30 . 
The panel (b) of the same figure shows that at t = ^2 , the 
coalescence instability among the secondary islands de- 
termines an inverse cascade till a three secondary-island 
system is reached. So far, the system still preserves its 
two-dimensional behavior. 

At this point, however, the system rapidly undergoes 
a full three-dimensional instability as shown in both pan- 
els of Fig. [3] referring to the instant t « 104 (ts). In the 
panel (a), the whole box is considered, while on the right 
a blow up is shown of the central region. The density 
structure is strongly kinked while the overall macroscopic 
configuration of the magnetic field topology is largely pre- 
served to its previous invariance along the Y direction. 
Several secondary instability modes gain ground onto the 
previously two-dimensional structure and drive a full space 
destabilization. The magnetic field loops around the high 
density regions following the plasma dynamics. The insta- 
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Fig. 4: Volume iso-surface of the density in the region 
of the current-sheet at the centre of our numerical do- 
main and magnetic field lines . The blue-red color ta- 
ble for the magnetic field lines describes the intensity of the 
X-component. The panel (a) shows the system at f ~ 108, 
whereas a detailed clip of the central island is shown in the 
panel (b). It is evident the final evolution (of the dominating 
mode) of the kinking instability driving the system and the 
complex configuration assumed by the magnetic field that is 
anchored at the structured. 

bility proceeds fast and in a few (dimensionless) crossing 
times of the current layer the central island is completely 
distorted by the kinking mode as shown in Fig. [4] panel 
(a). In the panel (b) of the same figure, a detailed view 
is presented of the central part of our numerical domain 
where the system is rearranged along the span-wise direc- 
tion (Y) with the magnetic field completely tangled to the 
plasma in a configuration highly unstable again to mag- 
netic reconnection. Indeed the system tears but in the 
orthogonal direction with respect the initial settings: re- 
connection spreads volume-wide. The system breaks and 
the disruption of the central island determines the expul- 
sion of higher-density regions transported by an increased 
velocity field as shown in Fig. [s] for t « 128 (t4). This 
fragmentation accelerates along the overall magnetic field 
lines and hits the two islands located at the side of the 
numerical domain. 

Mixing of the reconnected fluxes and fully tur- 
bulent state. — The turbulence and speed of the pro- 
cess of reconnection leads to an additional consequence: 
the destabilization of the downstream regions to a fur- 
ther Rayleigh- Taylor (RT) instability. As the plasma from 
the reconnection region is expelled, it collects downstream 
forming two regions of higher density, called plasmoids. As 
reconnection jets slow down and form the plasmoid, the 
acceleration field points contrary to the jet, while the jet 




Fig. 5: Volume iso-surface of the density in the region 
of the current-sheet at the centre of our numerical do- 
main, magnetic field lines, and velocity field. The sys- 
tem is considered at f ~ 112(^4), panel (a). The color tables 
for the magnetic field lines and the velocity arrow glyphs de- 
scribe the intensity of the x-component. The disruption of the 
central island produces the expulsion of higher-density regions 
transported by the increased velocity field. 

pushes towards the plasma that is piling up downstream. 
A configuration with acceleration and density gradient 
with opposite direction (the density increases along the 
reconnection jet, the acceleration points in the opposite 
direction) is unstable to the RT instability. The turbu- 
lence of the reconnection processes feeds this instability by 
providing the needed perturbations to allow the mode to 
grow. Conversely, the RT profoundly disturbs the 2D sym- 
metry facilitating the mixing in the 3D domain. Figure [6j 
panel (a)-(f), shows the characteristic fingers produced by 
the attempt of the heavy fluid to fall into the direction of 
the acceleration, just as a heavy fluid would do when lay- 
ered over a light fluid in a gravity field. The evolution of 
the instability completely destroys the overall structure of 
the system leading to a completely chaotic configuration 
as shown in Fig. [6] panel (f ) , where the density contour 
and magnetic field lines are shown at i « 160 (tg) • 
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